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Consideration of Moving Tooth 
Load in Gear Crack Propagation 
Predictions 

Robust gear designs consider not only crack initiation , but crack propagation trajectories 
for a fail-safe design . In actual gear operation , the magnitude as w ell as the position of 
the force changes as the gear rotates through the meslu A study to determine the effect of 
moving gear tooth load on crack propagation predictions was performed . Two- 
dimensional analysis of an involute spur gear and three-dimensional analysis of a spiral : 
bevel pinion gear using the finite element method and boundary element method were 
studied and compared to experiments. A modified theory for predicting gear crack propa- 
gation paths based on the criteria of Erdogan and Sih [18] was investigated. Crack 
simulation based on calculated stress intensity factors and mixed mode crack angle pre- 
diction techniques using a simple static analysis in which the tooth load was located at 
the highest point of single tooth contact was validated. For three-dimensional analysis , 
however, the analysis was valid only as long as the crack did not approach the contact 
region on the tooth. [DOI: 10.1115/1.1338118] 


Introduction 

Effective gear designs balance strength, durability, reliability, 
size, weight, and cost. Even effective designs, however, can have 
vthe possibility of gear cracks due to fatigue. In addition, truly 
robust designs consider not only crack initiation, but crack propa- 
gation trajectories. As an example, crack trajectories that propa- 
gate through the gear tooth are the preferred mode of failure com- 
pared to propagation through the gear rim. Rim failures will lead 
to catastrophic events and should be avoided. Analysis tools that 
predict crack propagation paths can be a valuable aid to the de- 
signer to prevent such catastrophic failures. 

Pertaining to crack analysis, linear elastic fracture mechanics 
applied to gear teeth has become increasingly popular. The stress 
intensity factors are the key parameters to estimate the character- 
istics of a crack. Analytical methods using weight function tech- 
niques to estimate gear tooth stress intensity factors have been 
developed [1,2]. Numerical techniques such as the boundary ele- 
ment method and finite element method have also been studied 
[3,4]. Based on stress intensity factors, fatigue crack growth and 
gear life predictions have been investigated [5-8], In addition, 
gear crack trajectory predictions have been addressed in a few 
studies [9-13]. 

From publications on gear crack trajectory predictions, the ana- 
lytical methods have been numerical (finite or boundary element 
method) while solving a static stress problem. In actual gear ap- 
plications, however, the load moves along the tooth, changing in 
both magnitude and position. No work has been done investigat- 
ing the effect of this moving load on crack trajectories. 

The objective of the current work is to study the effect of mov- 
ing gear tooth load on crack propagation predictions. Two- 
dimensional analysis of an involute spur gear using the finite el- 
ement method is discussed. Also, three-dimensional analysis of a 
spiral-bevel pinion gear using the boundary element method is 
discussed. A quasi-static numerical simulation method is pre- 
sented in which the gear tooth engagement is broken down into 
multiple load steps and analyzed separately. Methods to analyze 
the steps are discussed, and predicted crack shapes are compared 
to experimental results. 
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Two-Dimensional Analysis 

Gear Modeling. The two-dimensional analysis was per- 
formed using the FRANC (FRacture ANalysis Code) computer 
program developed by Wawrzynek [14]. The program is a 
general-purpose finite element code for the static analysis of two- 
dimensional cracked structures. The program uses principles of 
linear elastic fracture mechanics and is capable of analyzing plane 
strain, plane stress, or axisymmetric problems. A unique feature of 
the program is the ability to model cracks and crack propagation 
in a structure. A rosette of quarter-point, six-node, triangular ele- 
ments is used around the crack tip to model the inverse square- 
root stress singularity. Mode I and mode II stress intensity factors, 
K I and Kn , respectively, can be calculated using a variety of 
methods. (As a refresher, mode I loading refers to loads applied 
normal to the crack plane and tends to open the crack. Mode II 
refers to in-plane shear loading.) The stress intensity factors quanA 
tify the state of stress in the region near the crack tip. In the 
program, the stress intensity factors can be used to predict the 
crack propagation trajectory angles, again using a variety of meth- 
ods. In addition, the program has a unique re-meshing scheme to 
allow automated processing of the crack simulation. 

A spur gear from a fatigue test apparatus was modeled to dem- 
onstrate the two-dimensional analysis. The modeled gear had 28 
teeth, a 2ir pressure angle, a module of 3.175 mm (diametral 
piieii ui o in -1 ), and a face width of 6.35 mm (0.25 in). The gear 
had a backup ratio (defined as the rim thickness divided by the 
tooth height) of 3.3. The complete gear was modeled using mostly 
8-node, plane stress, quadrilateral finite elements. For improved 
accuracy, the mesh was refined on one of the teeth in which a 
crack was inserted. The total model had 2353 elements and 7295 
nodes. Four hub nodes at the gear inner diameter were fixed to 
ground for boundary conditions. The material used was steel. 

Tooth Loading Scheme. To determine the effect of gear 
tooth moving load on crack propagation, the analysis was broken 
down into eighteen separate load cases (Fig. 1). An initial crack of 
0.26 mm (0.010 in) in length was placed in the fillet of tooth 2, 
normal to the surface, at the location of the maximum tensile 
stress (uncracked condition). Six load cases were analyzed sepa- 
rately with the load on the tooth ahead of the cracked tooth, six on 
the cracked tooth, and six on the tooth after. The calculated stress 
intensity factors for unit loads at each of the load positions is 
shown in Fig. 2. These stress intensity factors were calculated 
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Fig. 1 Location of load cases for finite element mesh 


using the J-mtegral technique [15]. Loads on tooth 2 (cracked 
ooth) produced tension at the crack tip. The AT/s increased as the 
load moved toward the tooth tip (load cases 12 to 7, Fig. 2 b) due 
to the increased load lever arm. Loads on tooth 3 also produced 
tension at the crack tip; but at an order of magnitude less than 
those produced from the loads on tooth 2 (Fig. 2c). Loads on tooth 

(Ffg V 2n)° mPreSSi ° n t0 ^ CraCk aS shown ne 8 ativ e K/s 

actua * Ioacl magnitudes on the gear tooth were con- 
sidered as it went through the mesh. Computer program DANST 
(Dynamic ANalysis of Spur gear Transmission, [16]) was used for 
the analysis. This program is based on a four-degree-of-freedom 
torsional Jumped mass model of a gear transmission. The model 
includesdnving and driven gears, connecting shafts, a motor, and 
j load. The equations of motion for this model were derived from 
v^basic gear geometry, elementary vibration principles, and time- 
varying tooth stiffnesses. For simplicity, the static gear tooth loads 
of the solution were determined (Fig. 3). These loads were deter- 
mined from well-established gear tooth stiffness principles and 
static equilibrium. The loads are shown as a function of gear 
rotation for a driver torque of 68 N-m (599 in-lb). Tooth 2 began 
contact at a gear rotation of 10°. As the gear rotation increased, 
the load on tooth 2 gradually increased. Tooth 1 and 2 shared the 



Fig. 3 DANST computer program output of static gear tooth 
load, 68 N-m driver torque 


load for a rotation from 10 to 18°. From 1 8 to 23°, tooth 2 carried 
the complete load. At 23°, tooth 2 is considered at its highest point 
of single tooth contact (HPSTC). 

The stress intensity factors as a function of gear rotation were 
then determined by multiplying the stress intensity factors deter- 
mined from the units loads (Fig. 2) by the actual tooth loads (Fig. 
3) and applying superposition since linear elastic fracture mechan- 
ics was used. The results are shown in Fig. 4. As expected, the 
mode I stress intensity factor (Fig. 4a) was mostly influenced by 
the load on tooth 2. Note that the largest value of K t occurred at 
the HPSTC. Also note that the magnitude of AT; (Fig. 4a) was 
much larger than that of K n (Fig. 4b). This implied that K, was 
the driving force in the crack propagation, K n , however, affected 
the crack propagation angle as will be shown in the next section. 

Crack Propagation Simulation. From Williams [17], the 
tangential stress near a crack tip, cr^, is given by 
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Fig 2 Mode I and mode II stress Intensity factors for a ur 
load and an initial crack of 0.26 mm 


where r and ft are polar coordinates with the origin at the crack 
tip. Erdogan and Sih [18] postulated that crack extension starts at 
the crack tip and grows in direction of the greatest tangential 
stress. The direction of the greatest tangential stress is determined 
by taking the derivative of Eq. 1 with respect to ft setting the 
expression equal to 0 for, then solving ft Performing the math, 
this predicted crack propagation angle, 0 m , is given by 


ft m = 2 tan~ ! 



( 2 ) 


From Eq. 2, the predicted crack propagation angle is a function 
of the ratio of K l to K u . Erdogan and Sih [18] used brittle plexi- 
glass plates under static loading to validate their proposed theo- 
rems (i.e., the ratio of K f to Kfj was constant). For the gear prob- 
lem in the current study, however, the ratio of Kj to was not 
constant during gear rotation. This is shown in Fig. 4c (actually 
plotted as the ratio Kn to K t for clarity). In addition, Fig. 4 d gives 
the calculated 6 m from Eq. 2 as a function of gear rotation. 

In order to simulate gear crack propagation, a modification to 
the Erdogan and Sih theory was postulated in the current study. 
This modified theory states that the crack extension starts at the 
crack tip and grows in the direction of the greatest tangential 
stress as seen during engagement of the gear teeth. The procedure 
to calculate the crack direction is as follows: 1 ) K f and are 
determined as a function of gear rotation (Figs. 4 a and 4b , as 
described in the previous section), 2) the ratio of K t to K l{ as a 
function of gear rotation is determined (Fig. 4c), 3) 6 m (using Eq. 
2) as a function of gear rotation is determined (Fig. 4 d), 4) ctq 9 
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Ffg. 4 Stress intensity and tangential stress factors as a func- 
tion of gear rotation, 68 N-m driver torque, 0.26 mm initial crack 
size 


(using Eq. 1) as a function of gear rotation is determined (Fig. 4e), 
5) the predicte.d crack direction is the value of 6 m for which cr e$ is 
greatest during gear rotation. 

For t he gea r example given, the tangential stress factor (defined 
as cr qq\[ 2. Trr) is plotted as a function of gear rotation in Fig. 4e. 
This plot looks very similar to the mode I stress intensity factor 
plot (Fig. 4 a) since K t was much larger than K u (see Eq. 1). The 
tangential stress was largest at the HPSTC (gear rotation of 23°) 
and the predicted crack propagation angle at this gear rotation was 

Using this propagation angle, the crack was extended by 0.26 
mm (0.010 in), re-meshed, re-analyzed, and a new propagation 
angle was calculated using the method described above. This pro- 
cedure was repeated a number of times to produce a total crack 
length of 2.38 mm (0.094 in). The 0.26-mm crack extension 
length was based on prior experience in order to produce a smooth 
crack path. Figure 5 shows the stress intensity factors versus gear 
rotation for a number of crack lengths. Note that the mode I stress 
intensity factors looked similar but with increased magnitude as 
the crack length increased. In all cases, the selected crack propa- 
gation angles occurred when the tooth load was placed at the 




Fig. 5 Stress intensity factors from gear tooth crack propaga- 
tion simulation, backup ratio** 3.3 


HPSTC. Figure 6 shows a similar analysis but with a model of a 
thin-rimmed gear. Here, the gear was modeled based on the pre- 
vious design, but with slots incorporated in the rim to simulate a 
thin-rimmed gear. The backup ratio for this model was 0.2. As 
seen, the magnitudes of the mode I stress intensity factors during 
tension (gear rotations 18 to 45°) were larger than that of the 3.3 
backup ratio gear. Also, there was a significant increase in the 
compressive AT/s (gear rotation less than 18°) due to the increased 
compliance of the thin rim gear. 

Comparison to Experiments. Figure 7 shows the results of 
the analysis compared to experimental tests in a gear fatigue ap- 
paratus. The original model (backup ratio of 3.3) as described 
before was compared along with models of backup ratio of 1.0 




Fig. 6 Stress Intensity factors from gear tooth crack propaga- 
tion simulation, backup ratio —0.2 
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c) Backup ratio - 0.3. 


FI ?L 7 £° m P arlson Of predicted gear tooth crack propagation 
paths with experimental results (Predicted, E«experiments) 


and 0.3. These later two models were created using slots in the 
gear blank as previously described. The experiments were first 
reported by Lewicki and Ballarini [9]. Here, notches were fabri- 
cated in the tooth fillet region to initiate tooth cracking of test 
v gears of various rim thicknesses. The gears were run at 10,000 
rpm and at a variety of increasing loads until tooth or rim fracture 
occuiTed. As seen from the figure, good correlation of the pre- 
dicted crack trajectories to experimental results was achieved. For 
backup ratios of 3.3 and 1.0, tooth fractures occurred. For the 
backup ratio of 0.3, rim fracture occurred. 

As a final note, the analysis indicated that the maximum tan- 
gential stress at the crack tip always occurred when the tooth load 
was positioned at the HPSTC. Thus, for two-dimensional analysis, 
crack simulation based on calculated stress intensity factors and 
mixed mode crack angle prediction techniques can use a simple 
static analysis in which the tooth load is located at the HPSTC. 
This was based on a modification to the Erdogan and Sih crack 
extension theory and the fact that the mode I stress intensity factor 
was much larger than the mode II factor. 



Fig. 8 Boundary element model of OH-58 spiral-bevel pinion 


exported back to the FRANC3D program for post processing 
Fracture analysis, such as stress intensity factor calculation s, can 
then be performed. 

The spiral-bevel pinion of the OH-58C helicopter main rotor 
transmission was modeled to demonstrate the three-dimensional 
analysis. The pinion had 19 teeth, a 20° pressure angle, a 30° 
mean spiral angle, a module of 3.66 mm (diametral pitch of 6.94 
in ), and a face width of 31.75 mm (1.28 in). For OH-58 opera- 
tion, the pinion mates with a 71 -tooth spiral-bevel gear, operates 
at 6060 rpm, and has a design torque of 350 N-m (3099 in-lb). 

The boundary element model of the OH-58 pinion developed 
by Spievak [19] was used for the study. Three teeth, the rim cone, 
and the bearing support shafts were modeled (Fig. 8). The tooth 
surface and fillet coordinates were determined from the methods 
developed by Handschuh and Litvin [20] and Litvin and Zhang 
[21]. The mesh of the three teeth was refined for improved accu- 
racy. A half-ellipse initial crack with major and minor diameters 
of 3.175 and 2.540 mm, respectively (0.125 and 0.100 in), was 
placed in the fillet of the middle tooth normal to the surface. The 
crack was centered along the face width and centered along the 
fillet. The complete gear model had a total of about 2600 linear 
elements (both triangular and quadrilateral) and about 2240 nodes. 
For boundary conditions, the end nodes of the larger diameter 
shaft were fixed and the nodes on the outer diameter of the smaller 
diameter shaft were constrained in the radial directions. Again, the 
material was steel. 


Three-Dimensional Analysis 

Gear Modeling. The three-dimensional analysis was per- 
formed using the FRANC3D (FRacture ANalysis Code for 3 Di- 
mensions) computer program developed by Wawrzynek [14]. This 
program uses boundary element modeling and principles of linear 
elastic fracture mechanics to analyze cracked structures. The ge- 
ometry of three-dimensional structures with non-planar, arbitrary 
shaped cracks can be modeled. The modeling of a three- 
dimensional cracked structure is performed through a series of 
programs. Structure geometry grid point data are imported to a 
solid modeler program. Here, appropriate curves and faces (or 
patches) are created from the grid data as well as a closed-loop 
s “ rl erf geometr y modeI * This surface model is then imported to 
the FRANC3D program for boundary element model preparation. 
The user can then mesh the geometry model using 3 or 6 node 
triangular surface elements, or 4 or 8 node quadrilateral elements 
Boundary conditions (applied tractions and prescribed displace- 
ments) are applied on the model geometry over faces, edges, or 
points. Initial cracks, such as elliptical or penny shaped, can be 
inserted in the structure. After complete formulation, the model is 
shipped to a boundary element equation solver program. Once the 
displacement and traction unknowns are solved, the results are 


Tooth Contact Analysis and Loading Scheme. Due to the 
geometrical complexities and three-dimensional action, numerical 
methods are required to determine the contact loads and positions 
on spiral-bevel teeth since no closed-form solution exists. The 
method of Litvin and Zhang [21] was used to determine the mean 
contact points on the spiral-bevel pinion tooth. The method mod- 
eled tooth generation and tooth contact simulation of the pinion 
and gear. With the mean contact points taken as the centers, con- 
tact ellipses were determined using Hertzian theory [22]. Figure 9 
shows the estimated contact ellipses on the spiral-bevel pinion 
tooth. Fifteen separate ellipses (load cases) were determined, start- 
ing from the root of the pinion and moving toward the tooth tip 
and toe. Load cases 1-4 and 12-15 were double tooth contact 
regions while load cases 5-11 were single tooth contact regions. 
Note that load case 1 1 corresponds to the load at the HPSTC. For 
each load case using the boundary element method, tractions were 
applied normal to the surface to the appropriate ellipse with the 
magnitude equal to the tooth normal force divided by the ellipse 
area. 

Crack Propagation Simulation. The procedure for the three 
dimensional crack propagation simulation of the OH-58 spiral- 
bevel pinion was as follows. For each of the load cases of Fig. 9, 
the mode I and mode II stress intensity factors were determined at 
25 points along the crack front (note that for three-dimensions, 
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Table 1 Results of multiple load case crack simulation 
analysis. 



Fig. 9 Location of tooth contact ellipses and magnitude of 
load on OH-58 spiral-bevel pinion tooth 


there is a crack front, not just a crack tip as in two-dimensions). 
v The extended crack directions at each of these 25 points were 
determined using the modified Erdogan and Sih crack extension 
theory as described in the two-dimensional analysis. That is, as 
the cracked spiral-bevel pinion tooth was engaged in the mesh, the 
crack extension started at each point along the crack front and 
grew in direction of the greatest tangential stress at those points 
during mesh. The amount of crack extension at each point along 
the crack front was determined based on the Paris crack growth 
relationship [23] where 


KlJ 


K 


/.max, 


(3) 


where a , was the amount of extension of the i lh point along the 
crack front, Kj , was the mode I stress intensity factor of the i th 
point along the crack front corresponding to the load case which 
gave the largest tangential stress for that front point, K } max was 
the value of the largest K f i along the crack front, a max was the 
maximum defined crack extension along the crack front, and n 
was the Paris material exponent. From experience, the maximum 
extension size, a tMXt was set to 1.27 mm (0.050 in.). The Paris 
exponent, /i, was set to 2.954 based on material tests for AISI 
9310 steel by Au and Ke [24]. A third-order polynomial was then 
used to smooth the extended crack front. The new crack geometry 
was then re-meshed. After remeshing, the model was rerun and 
solved for stress intensity factors, and crack propagation direc- 
tions. The above procedure was repeated a number of times to 
simulate crack growth in the gear tooth. 

Table I gives results from the first four steps during this pro- 
cess. Note that step 0 corresponds to the initial half-ellipse crack. 
For steps 0 and 2, the largest tangential stress occurred at the 
HPSTC (load case 11) for the majority of the points along the 
crack front. For steps 1 and 3, the largest tangential stress oc- 
curred at load cases 8, 9, 10, or 1 1. 

Figure 10 shows the stress intensity factor distribution along the 
crack front for step 1 (crack area of 5.96 mm 2 (0.009 in 2 )). Similar 
to the spur gear analyses, the K{ s were larger as the load moved 
from the root to the tip due to the larger load lever arm. Other than 
absolute magnitude, the K t distributions along the crack front 
looked similar for the various load cases. Figure 1 1 depicts the 


Step 

Crack area 
(mm 2 ) 

Crack front point(s) 

Load case for 
largest a 89 

0 

3.12 

1 

8 



2-25 

11 

1 

5.96 

1 

9 



2,4-7,21,23-25 

10 



3,8-20,22 

11 

2 

10.35 

1 

10 



2-25 

11 

3 

13.35 

7,9,20 

8 



5-7,10-15,21,26,27 

9 



1-4,8,16-19,22-25 

11 


stress intensity factors plotted against load case (at a point along 
the front, biased toward the toe, normalized position along the 
crack front of 0.83). This figure shows the simulated distribution 
as the pinion engages in mesh with the gear. Note again that the 
ratio of K u to Kj was not constant during engagement. 

Figure 12 shows exploded views of the pinion crack simulation 
after seven steps. It should be noted that the loading was placed 
only at the HPSTC for the last three steps. This was due to mod- 
eling difficulties encountered using the multi-load analysis. It was 
felt that this simplification did not significantly affect the results 
due to the smoothing curve-fit used. In addition, the tangential 
stress near the crack tip was either largest, or near its largest 
value, when the load was placed at the HPSTC. 

Comparison to Experiments. Figure 13 shows the results of 
the analysis compared to experimental tests. The experimental 
tests were performed in an actual helicopter transmission test fa- 
cility. As was done with the gear fatigue tests described before, 
notches were fabricated in the fillet of the OH-58 pinion to pro- 




Normalized position along crack front 


Fig. 10 Stress intensity factors from three dimensional OH-58 
pinion tooth crack propagation simulation; step 1, crack area 
-5.96 mm 2 
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Fig. 1 1 Stress intensity factors from three-dimensional OH-58 
pinion tooth crack propagation simulation; step 1, normalized 
position along crack front- 0.83 


y 

mote fatigue cracking. The pinion was run at full speed and a 
variety of increasing loads until failure occurred. Shown in the 
figure are three teeth that fractured from the pinion during the 
tests (Fig. 13 b). Although the notches were slighdy different in 
size, the fractured teeth had basically the same shape. 

A side view of the crack propagation simulation is shown in 
Fig. 1 3 a for comparison to the photograph of the tested pinion in 
Fig. 13 b. From the simulation, the crack immediately tapered up 
toward the tooth tip at the heel end. This trend matched that seen 



Fig, 12 OH-58 spiral-bevel pinion tooth crack propagation 
simulation after seven steps 



7 a) Simulation. ' \ V ^ : 1:’ b) E>£er]^ 


Fig. 13 Comparison of OH-58 spiral-bevel pinion tooth crack 
propagation simulation to experiments 


from the tests. At the toe end, the simulation showed the crack 
progressing in a relatively straight path. This also matched the 
trend from the tests. Toward the latter stages of the simulation, 
however, the crack tended to taper toward the tooth tip at the toe 
end. This did not match the tests. One problem encountered in the 
simulation during the later steps was that the crack at the heel end 
of the tooth became close to the actual contact ellipses. It was felt 
that the crack-contact interaction may have influenced the trajec- 
tory predictions to cause the discrepancy. 

Spievak [19] reported on another method to account for the 
non-uniform K n to Kj ratio during pinion tooth engagement. This 
method considered contributions from all load cases in the crack 
angle prediction scheme and presented a method to cumulate the 
load effects. From these studies, reported crack propagation simu- 
lation of an OH-58 pinion also predicted the erroneous taper to- 
ward the tooth tip at the toe end. Again, the crack-contact inter- 
action may have influenced the trajectory predictions to cause the 
discrepancy. Spievak also reported on a simulation using only the 
load at the HPSTC. The results from that simulation were similar 
to that achieved in the current study. It should be noted that the 
proposed method in the current study to account for moving tooth 
loads for the three-dimensional analysis was extremely cumber- 
some. It is therefore felt that the analysis using only the load at the 
HPSTC appeared accurate as long as the crack did not approach 
the contact region on the tooth. 


Conclusions 

A study to determine the effect of moving gear tooth load On 
crack propagation predictions was performed. Two-dimensional 
analysis of an involute spur gear using the finite element method 
was investigated. Also, three-dimensional analysis of a spiral- 
bevel pinion gear using the boundary element method was dis- 
cussed. The following conclusions were derived: 

1) A modified theory for predicting gear crack propagation 
paths based on the criteria of Erdogan and Sih was validated. This 
theory stated that as a cracked gear tooth was engaged in mesh, 
the crack extension started at the crack tip and grew in direction of 
the greatest tangential stress during mesh. 

2) For two-dimensional analysis, crack simulation based on cal- 
culated stress intensity factors and mixed mode crack angle pre- 
diction techniques can use a simple static analysis in which the 
tooth load is located at the highest point of single tooth contact. 

3) For three-dimensional analysis, crack simulation can also use 
a simple static analysis in which the tooth load is located at the 
highest point of single tooth contact as long as the crack does not 
approach the contact region on the tooth. 
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